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Abstract 

We present an iterative algorithm to compute numerical approxima¬ 
tions of the potential for the Schrodinger operator from scattering data. 
Four different types of scattering data are used as follows: fixed energy, 
fixed incident angle, backscattering and full data. In the case of fixed 
energy, the algorithm coincides basically with the one recently introduced 
by Novikov in [18], where some estimates are obtained for large energy 
scattering data. The numerical results that we present here are consistent 
with these estimates. 


1 Introduction 


We consider the following scattering problem appearing in quantum physics, 
f (—A + V{x) — k^)u{x) = 0, a; G R^, 

\ U = Ui + Ua ^ 


where Ui = Ui{x,9,k) = ® is the incident wave, with wave number 

A: € R and direction of propagation 9 G S^. The potential V {x) € Z/°°(R^) 
is assumed to be compactly supported and Ua{x, 9, k) is the scattered wave 
that satisfies the Sommerfeld radiation condition at infinity, 

‘^^{x,9,k) — ikua(x,9,k) — as r —>■ oo, (2) 

or 


with r = \x\. From this condition we can deduce the following asymptotics 

for Ua, 


Ua{x,9,k) 


.if,e,k) 




1 / -l/2\ 

-l-o(r ' ) 


( 3 ) 


where the function Uaa{9', 9, k) depending on fc G R, the incident angle 9 G 
and the reflecting angle 9' £ S^, is known as the scattering amplitude 
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or far field pattern and constitutes the data of the inverse problem (see 
Ch. 5 in [8] for details). 

We are interested in recovering the quantum mechanical potential V [x) 
from the far field data Uoo{0', 9, k) for some values ot 9' & , 9 & and 

k gR. 

Let us state the problem in an equivalent integral formulation. Under 
the above conditions on V{x), the system (l)-(2) is equivalent to the so- 
called Lippmann-Schwinger equation 

u{x,9,k) = + [ ^{k\x - y\)V{y)u{y,9,k) dy (4) 

where <l?(r) = \H^\r) and is the Hankel function of the first kind 
of order zero. On the other hand, the potential and the far field pattern 
are related by the integral equation (see [9], [8]), 

Uoa{d',9,k)= [ V{x) u(x,9,k) dx. (5) 

Therefore, the inverse scattering problem can be stated as follows: Given 
the far held data Ucx>i9',9,k), for some values of 9' G 9 G and 
k G M, hnd a compactly supported V(x) G that satishes (5), 

where u is the solution of (4). 

Note that Uoo depends on three variables while V only depends on two. 
Thus, in principle, V(x) could be recovered from the partial knowledge of 
Uoa{9',9,k) on a suitable 2-dimensional submanifold of x 5^ x R. The 
usual choices are: hxed energy k = ko, hxed incident angle 9 = 9o and 
backscattering 9' = —9. 

A number of algorithms have been proposed to recover the potential 
from the scattering data. Most of them are based on the Neumann-Born 
series for V described in [12] (see also [23] and [24]) which roughly consists 
in substituting (4) into the right hand side of (5) iteratively. At each step a 
new multilinear term in V appears in the right hand side of (5). Assuming 
that V can be written in power series with respect to a small parameter 
e, i.e. V = making equal the same powers of e one 

easily obtains an iterative formula for U". The first term of the series 
provides a somehow linearization of the inverse problem and it is known 
as the Born approximation of the potential. Other strategies are based 
on perturbation methods (see [6]) or more direct inversion algorithms as 
those proposed by Novikov ([14], [15], [16], ]17] and [19]). These have been 
implemented numerically in [3] and [5], (see also the references therein). 

In this paper we investigate numerically an iterative algorithm to ap¬ 
proximate the potential based on a suitable hxed point iteration on the 
integral formula that dehnes the Born approximation, which we denote 
by Vb- We obtain a sequence of approximations where the new 

approximation can be deduced from U" and its associated scattered 

held u" by solving a single inverse Fourier transform. To validate the al¬ 
gorithm we consider a numerical approximation of the inversion formula 
for the Fourier transform. A convergence result for this numerical method 
to approximate the inverse Fourier transform is given in Section 4, using 
a technique which may be understood as the Nyquist-Shannon sampling 
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theorem applied to the Fourier transform. As applications we obtain, on 
one hand, estimates for the numerical approximation of the Born approx¬ 
imation and, on the other hand, numerical evidences that the sequence 
converges and it provides a good approximation of the potential 
V(x) in few iterations. 

For scattering data at fixed energy k = ko, the algorithm is similar 
to the one introduced by Novikov in [18], where the convergence of the 
iterative process is investigated from the theoretical point of view. The 
main difference with the algorithm described here is that in [18], at each 
iteration, the approximation of the potential is modified by a low-pass 
filtering/cutting process that we briefly explain at the end of Section 2 
below. The convergence result in [18] is stated for any dimension d > 2 
and for smooth potentials V (more than d derivatives). The numerical 
results described here suggest that the convergence result in [18] could be 
true also for less regular potentials {V G L°°) and without the filtering¬ 
cutting process. 

The rest of the article is organized as follows. The iterative algorithm 
to recover the potential is given in Section 2. The numerical method to 
approximate both the Born approximation and the subsequent approx¬ 
imations V" is described in Section 3. An analysis of the convergence 
for the numerical approximation of Vb in terms of the mesh step is per¬ 
formed in Section 4. In Section 5 we show how the scattering data was 
simulated from a potential example, and numerical experiments for both 
the Born approximation and the sequence 14”. Finally, some conclusions 
are presented in Section 6. 


2 The iterative algorithm to approximate 
the potential 


In this section we present the iterative algorithm to approximate the po¬ 
tential. To introduce it we first describe the Born approximation in detail. 
When substituting the equality u = Ui-\-Us into the right hand side of (5) 
we obtain. 




+ f e V{x) Ua{x,6,k) dx. (6) 

The Born approximation Vb to V is dehned formally as the solution of 
(6) when neglecting the last nonlinear term, i.e. the solution of the linear 
problem. 



ilR2 

However, this identity is not consistent. Given ^ G R^, the right hand 
side is constant for those values (6', 6, k) in the set 

Gj = {{O', e, k)e X X R, such that k{e' -9) = 2tt£}, (8) 
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while the left hand side of (7) does not satisfy this compatibility condition 
necessarily. Therefore, a proper definition of the Born approximation Vb 
requires also a strategy to select, for a given ^ £ R^, some specific values 
{d '£ Gg in a unique way. 

Following the methodology in reference [22], in order to define Vb one 
has to choose two open subsets M, of 5"^ x 5^ x R and R^, respectively, 
such that the operator 0 : M i—> SI defined by <^(0^ 0, k) = (l/27r)fc(0' 
is an isomorphism. This way, from the definition (7), it follows that 
Ucx>{4>~^(0) — for C £ where T denotes the Fourier transform 

operator. Hence, provided that R^ \ SI has zero measure, one can easily 
obtain Vb by inverting the Fourier transform, 

Vb{x) = [ 

Jn 

= f (9) 

Jn 


Given ^ = (^ 1 ,^ 2 ) £ SI, we determine the values ^(?)i ^(0) ^ 

for the usual types of inverse potential scattering: 

1. Case 1: Fixed energy Vb- We fix fc > 0. Here, SI = D(0, fc/7r), 
M = {{9',e,k) ■ 9,9' G 5^}. Define (9'{^),9{^),k) £ as follows. 


0(C) = 0(Q!)||j, where <d{a) 

9'{0^ei0 + 2nl 


COS a — sin a 
sin a cos a 


( 10 ) 


Here, ©(a) performs the rotation of angle ce{^) = arccos(—rr^) with 
respect to the origin. Note that, in this case, we can only describe 
values of C in the disk |C| < fe/yr. Indeed, the set R^ \ S2 does not 
have zero measure. Thus, in order to define properly the Fourier 
transform of Vb, we extend the far field pattern Uaa{9', 9, k) by zero 
outside M = ())“^(S1), which means that a low pass filter is applied 
to the Born approximation. In this case, Vb can not be compactly 
supported. 

2. Case 2: Fixed incident direction Vb- We fix 0 £ 5^, so H = {C € 
R^ : C'0 7^ 0} andM = {(0',0,A:) : 0' £ S^\{0}, k £ R\0}. Define 
fe(C) and 0^(C) as follows. 


3 . 


r fc(C = -^l^l7(C-^), 

10'(e) = 7c) + 2^4y. 

Note that, in this case, we can only describe values 
hyperplane C • 0 = 0, but this is a zero measure set. 

Case 3: backscattering Vj^“® . We fix 9' = —0, so 
choose M = {(—0,0, A:) : 0 £ S^, fc > 0} (the case 
possible as well). Define fc(C) and 0(C) as 

r fc(c) = ICK, 

I = 


(11) 

of C outside the 

H = R^ \ 0, and 
fc < 0 would be 

( 12 ) 
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4. Case 4: full scattering data Vg. In this case the Born approximation 
is defined as an average of the Born approximations for fixed incident 
angle, i.e. 

viix) = ^ de. (13) 

In general, all these Born approximations recover some properties of 
the potential V{x) and, in particular, its singularities in the scale of 
weighted L^, Holder and Sobolev spaces, see [22], [20], [26], [27], [25], 
[1] and [2]. 

Once described the Born approximation to V we introduce the itera¬ 
tive algorithm to improve this approximation, which is basically a fixed 
point iteration applied to (6). We consider the sequence of approximate 
potentials defined recursively by, 

' = 0 , 

For n = 0,1, 2,1/"+^ is the solution of 
' uoo{0',e,k) = J^,e-*'^^‘>'-<>>^V^+\x)dx 

. +J^2e-^'‘'>'-^V"ix)u:ix,e,k)dx, 

where u'^{x,0,k) is the scattered field associated to the potential V", i.e. 
u'^ = u" — ® and u" is the solution of (4) with V = V". To simplify, 

we derive directly w" as the solution of the following Lippmann-Schwinger 
equation, 

u:(x,e,k) = f ^k\x-y\)V-{y)e^^^Uy 

+ f ^{k\x-y\)V"{y)u"{y,e,k) dy, (15) 

Jr2 

which is equivalent to (4). 

Note that satisfies the equation (9) for the Born approximation 
(i.e. = Vb) whereas, for n > 1, V" can be interpreted as an improved 

approximation of V, since it takes into account a better approximation of 
the second integral in the right hand side of (6). Similar iterative proce¬ 
dures have been previously used by A. Ruiz in [26, Section 5] for theoret¬ 
ical issues, where an alternative second approximation is constructed 
by plugging the Born approximation into the first nonlinear term of the 
Neumann-Born series, or by Y.M. Chen and W.C. Chew in [7] where a 
similar idea is considered to improve the Born approximation in a related 
electromagnetic inverse scattering problem. 

As for the Born approximation, the system (14) is not consistent in 
general and a suitable strategy is required to recover from (14)-(15). 
To be more precise, the new potential 1/”+^ is obtained from the previous 
one 17" through the Fourier transform T, 

-[ ^^^"^V'^{y)ul{y,0{i),k{i)) dy, (16) 

where ^iO) ~ i® defined according to one of the strate¬ 

gies defined before (fixed k, fixed 9, etc.), and ^ G 17. In this way, we have 
a different sequence {17"}n>o for each one of the cases defined before. 


5 



The iterative algorithm introduced above is far from being justified 
from the theoretical point of view. A rigorous validation would require to 
address in particular the following issues: 

11. The existence of 

12. If II holds true, the convergence of the sequence {V^}n>i. 

13. If 12 holds true, the identification of the limit V* with the original 
potential V. 

Concerning the first issue, the definition of with n > 1, requires 

previously to compute the solution u" to the Lippmann-Schwinger equa¬ 
tion (15) corresponding to I/". We do not know whether V" satishes the 
conditions guaranteeing the existence and uniqueness of solution to this 
equation. 

The second issue is also completely open. The third one is related 
to the uniqueness for the potentials V(x) satisfying (6) but only when 
{9', 6, k) satisfies the restrictions provided by the case we are considering. 
For instance, in the case of the backscattering {O' = —9), assuming that 
we can pass to the limit in (14), V* will satisfy 

Uoo{6',9,k)= [ ''''V*{x) u*{x,9,k) dx, (17) 

7*2 

when 9 = —9' and A: > 0. Here u* is the solution of (4) with the potential 
V*. Therefore V{x) = V*{x) if equation (17) determines V uniquely. 
This is also a difficult question, although local and generic uniqueness 
have been proved by Eskin and Ralston ([10]). The uniqueness problem 
for fixed energy has been solved for Bukhgeim [4] (see also [11]) and for 
fixed incident direction is also an open problem (see [29] and [28] ). If we 
use all scattering data, Uoo{d', 9, k) determines the potential uniquely (see 
[28] and [21]). 

As we said in the introduction, the algorithm presented here is closely 
related to the one introduced in [18] for fixed energy scattering data. 
In this case, the sequence of approximations {IF"}„>o starts also from 
= 0 and it is constructed iteratively as follows: we consider again for¬ 
mula (16) but replacing V” by W" and by IF"+^, an intermediate 

function. The new iteration is then obtained from in two 

steps: we first compute a suitable low-pass filter of and then mul¬ 

tiply it by a compactly supported cutoff function. With this strategy the 
issue II above is satisfied for {IF"}„>o, due to the fact that IF” is com¬ 
pactly supported at each step (by construction), while the convergence of 
the sequence {IF”}„>o is deduced from a careful analysis based on the 
filtering step. In this way, the following estimate is obtained, 

||IF”-F[[l=o (18) 


where Cn is uniformly bounded. 





r — ci\”\ r — d 
r ) ) 2d 


d is the space dimension (in our case d = 2), and r > d is the required 
number of derivatives in of the potential V. Note that this estimate 
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does not imply the convergence of W”' to V as n —>■ oo, for fixed k, but 
it provides a polynomial convergence rate as A: —>■ oo, for fixed n. We 
compare below our numerical estimate for — V with this one. 


3 Numerical approximation 

In this section we describe the numerical approach to find approximations 
of both the Born approximation dehned by formula (9) and the sequence 
of iterative approximations defined by (16). We hrst introduce some 
notation on the hnite dimensional space and then we state the numerical 
versions of (9) and (16). 


3.1 Finite dimensional trigonometric space 

Given i? > 0, we dehne 

Gr = {x = {xi,X 2 ) € : \xk\ < R, A: = 1, 2} . 

The family of exponentials 

‘Piix) = 0 = ^, j = Uij2)eZ^, 

constitutes an orthonormal basis on L^{Gr) with the norm 


II l|2 

ll«llo = 


i2Ry 


Ij 


\u{x)\ dx. 


We also introduce the space H\ = H\(Gr) which consists of 2i?—multiperiodic 
functions (distributions) having finite norm 


1/2 


with 


and 


w U 




,iez2 




G 5 (0, 0) 7^ j € I? 

1, i = (0,0), 


M(^i) = / u{x)ifij{x)dx, 

Jgr 

the Fourier coefficients of u. 

We now introduce a hnite dimensional approximation of Hx. Let us 
consider h — 2R/N with N £ N and a mesh on Gr with grid points jh, 
j G Zl and 

Z" = {i = (ii,i2) G < ife < y, fc = 1,2| . 

We also consider ‘Jh the Hnite dimensional subspace of trigonometric poly¬ 
nomials of the form 


Vh= Cjifj, 
jezl 


Cj G C. 


7 




Any Vh € 7h can be represented either through the Fourier coefficients 


Vhix) = Vh{j) 

or the nodal values 

Vh{x) = Y "^hUh) ^h,j{x), 


where 

/ \ l2 X ' iTrk-ix—jh) / R 

Vh,j{x) = h 2 ^ e '■ J n , 
kezl 

For a given vh £ Th, the nodal values Vh and the Fourier coefficients Vh 
are related by the discrete Fourier transform Ih as follows, 

2 - - 1-1 
Vh = h JhVh, Vh = Vh, 

where, as usual, 3^h relates the sequence x{k) with X{k) according to 


iV-l 

X{k) = Y , fc = 0,1,..., N-1. 

n=0 

The orthogonal projection from Hx to 7h is defined by the formula 


PhV = Y 

while the interpolation projection QhV is defined, when A > 1, by 
QhV £ 7h, (Qhv)(jh) = v(jh), j € Zl. 


3.2 Finite dimensional setting 

Let R > 0 such that supp(y) C Gr. We approximate the Born approxi¬ 
mation in (9) by the following finite dimensional version: Find Vs,h G 7h 
such that 


h^7hVB,hU) = «oo(0'fe),e(G),fc(G)), j e (19) 

Note that Ve.h is computed from a single inverse discrete Fourier trans¬ 
form from the values of the far field Uoo{9'{-), d{-), k{-)) at the mesh points 
with j G More precisely, the numerical approximation is obtained 
from the following process: 

Algorithm 1: 

1. Choose h according to the mesh grid where we will compute the 
nodal values of VbX- Xj = jh, j G 

2. Construct the mesh = j/{2R) with j G 

3. Compute Uoa{d'{■),9{-),k{-)) at the mesh points ^j. 


4. Invert the discrete Fourier transform to obtain the values of Vbm at 
the nodes Xj. 

Analogously, we compute the numerical approximations of 1/"+^ in 
(16) by solving the following discrete formula for G 7h 

k&zl 

where u"{x,0{^j),k{^j)) is obtained from a discrete formulation of the 
Lippmann-Schwinger equation (15). In the experiments that we present in 
Section 5 below we follow the trigonometric collocation method introduced 
by G. Vainikko [30]. 

The iterative process is then as follows: 

Algortihm 2: 

1. Choose h according to the mesh grid where we will compute the 

nodal values of VbX- J ^ 

2. Construct the mesh = j/{2R) with j £ Z^. 

3. Compute Vb.h =’■ following the Algorithm 1 above. 

4. For n = 1,..., K-. 14’*+^ 

(a) Compute u"(x,9{^j),k(^j)) by approximating the Lippmann- 
Schwinger equation (15) with V = V))*, for each mesh point 

j ^ 

(b) Compute the right hand side of (20) for each mesh point ^j. 

(c) Invert the discrete Fourier transform to obtain the values of 

at the nodes Xj. 

5. End 

Note that each iteration requires N'^ performances of the Lippmann- 
Schwinger equation solver, which is a costly computation. The good news 
is that, on a hand, these computations can be parallelized, and on the 
other hand, the results can be estimated by interpolation with the ones 
generated from a few computations on a coarse grid. 

In the next section, we give estimates for the error of the numerical 
approximation VB,h in (19) with respect to a periodized version of the 
Born approximation Vb- As a result, if Vb has compact support in Gr, 
we have estimates for the difference between the approximation VB,h and 
the true Vb, since in this case the periodized version coincides with the 
Born approximation Vb in Gr. We do not know whether the Born ap¬ 
proximation preserves or not the compactness of the potential’s support, 
except for the fixed energy case, in which a low pass filter is applied to 
Vb, and consequently Vb is not compactly supported. 
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4 Convergence of the numerical approx¬ 
imation of the Born approximation 

In this section, we prove a convergence result for the numerical approxima¬ 
tion of the Born approximation obtained by the Algorithm 1. The result 
is a direct consequence of a Lemma below that provides an estimate for a 
discrete approximation of the inverse Fourier transform. 

Theorem 1. Let Vb the Born approximation of a potential V{x) defined 
by (9), according to one of the possible cases 1-4 defined in the introduction 
(fixed energy, fixed incident angle, etc.). Consider a mesh of size h as 
before in such a way that € fl for all ji € Let Vg the periodized 
version of Vb defined as 

E + (21) 

If Vg £ for some A > 0, then 

\\VB.H-Vl\\o<h^\\V^h, ( 22 ) 

where Vb.h £ “Lh is the solution of (19). 

Remark 2. The periodie version of Vb given in (21) will coincide with 
Vb in Gr only if Vb is compactly supported in Gr. We do not know if 
this is the case, in general. 

The proof of Theorem 1 is a direct consequence of the following lemma, 
which is an adapted version of the sampling theorem. 

Lemma 3. Let q(x) € L^(R^) be such that 

where g{fi) is a continuous known function. Consider the finite dimen¬ 
sional approximation given by qu £7h, the solution of 

h^7hqh{j) = g{(,j), j G Z|. (23) 

Let 5 ** be the periodized version of q given by 

g“(®) = E ^(a: + “^Rj)- (24) 


In particular, if q is compactly supported in Gr then q^ = q on Gr. 

Then, if g** £ for some A > 0, 

Ikh - g^llo < h^||(7*'||A- (25) 

Proof: We divide the proof in several steps. 

Step 1: Filtering. In this first step we dehne a low pass hlter of g{fi), 
denoted by gLp{^), in such a way that 

supp {'J~^gLp) C Gr. 
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We take 


where 


5LP(C) = - Cj), 

J6Z2 


sin(27r_R5^) sin(27r7?^^) 
2ttR^^ 2-kR^^ 


with ^ 


and = j/{2R), j £l^. Note that with this choice 
gLp{ij) = g{ij), j e 
Let us define q £ as the solution of 

ya(e) = flpp(o, 


{e,e), 


(26) 


By construction q has compact support in Gr. Moreover, taken into 
account that the Fourier transform of is (2_R)^m(f — ^j) 

we have, 


q{x) = 3" ^gLp{x) 


1 

{2Ry 




(27) 


where XG^ (a;) is the characteristic function of the two dimensional interval 
Gr. 


Step 2: Cut and periodize. The main idea is to establish a convergence 
result for q in Gr with the Hx metric. Therefore we define q^ £ Hx as the 
2i7-periodic extension of q which coincides with q in Gr. We show that 
Phq^ is in fact qh defined by (23). 

First of all, write the restriction of q{x) to Gr in Fourier series, i.e. 

q{x) = q\x) = ^ Cjvaj(x), x £ Gr. (28) 

jez2 

It is easy to see, from equation (27), that Cj = gLp{^j) = g(Cj)- 
The projection Phq^ £ '^h is given by 

Phq\x) = ^ Cjipj{x), a: G K^. 

aez? 

Then, we have for the nodal values of Pnq, Phq 

h^^h{lW)U) = g{Q, j&^l, 

which is the formula that defines qh- Therefore, 

qh = Phq^. 


Step 3. Now we estimate the norm in (25) 


gh - = \\Phq^ - q^llo < {4:h)^\\q^\\x, 
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where we have used that the projection operator Ph satisfies, in general, 


\\PhV - nil A < (4h)'^ 


for all V G fJ. > X- 


Step 4. Finally, we prove that q** is given by formula (24). As q = 
in Gr it is sufficient to prove it for q. Taking into account (27) we have, 



<l{x + ‘2Rj)xGn{x) 


i6Z2 


where we have used the Poisson summation formula in the last identity. 
This completes the proof of (24). 

5 Numerical experiments 

In this section we show the efficiency of Algorithms 1-2 to approximate 
the potential, when inverting the Fourier transform with the numerical 
method described in the previous section. We divide this section into 
three subsections. The first one is devoted to the direct problem, i.e. we 
show how to construct synthetic far field data for a given potential in 
order to test our algorithms. This is also useful to apply the Algorithm 
2. The approximation of the potential from the far field data is discussed 
in the next two subsections, we consider separately the results with the 
Algorithm 1, that provides the Born approximation, and the Algorithm 
2 . 

5.1 Data simulation 

We consider different situations according to the different cases presented 
in the introduction. The far field data Uoo (0^(C)i ^(C)> is simulated by 
solving the direct problem with the numerical code written by K. Knudsen, 
J. L. Mueller and S. Siltanen to solve numerically the D-bar equation for 
the two-dimensional Calderon problem. These codes are based on the 
numerical approach to solve Lippmann-Schwinger equations introduced 
by G. Vainikko [30]. 

For each scattering type, consider the mesh {^j = j/{2R) : j € Z\}. 
Fix j and write Oj := kj := k{^j). We compute a numerical approx¬ 

imation Us{xi,9j,kj) to the solution Us{-,6,k) of the problem (l)-(2) for 
9 = 9j, k = kj, hy adapting the aforementioned codes by S. Siltanen et 
al. to a 2i?-biperiodic version of the integral equation (4). To this end, 
we write (4) as follows 


[7 - 4-fc * (F • (•))]^^s(-, 0, k)^$k* 


(29) 


where 4>fc(x) := <I>(fc| 2 ;|) denotes the fundamental solution to the Helmholtz 
equation, and periodize the equation (29) by replacing all its terms with 
their 27?-biperiodic extension, cutting the Green’s function smoothly and 
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taking the convolution operator on the torus. We follow the spirit of 
Section 15.4 in the book [13] by J. L. Mueller and S. Siltanen. The funda¬ 
mental solution i&fc is computationally implemented through the Matlab 
function besselh. 

The approximation Us{xi, 9j, kj) is computed on a grid {xi — Ih : I £ 
Z|} of the square Gr, with h = 2R/N. Assuming a potential supported in 
the unit disc D{0, 1), the condition i? > 2 is necessary for the periodization 
argument. 

We have added a 5% Gaussian noise to the synthetic scattering data 
Ug{xi,9j,kj) to simulate possible measurement errors and validate the 
robustness of the approach. More precisely, we take (1 -I- 0.05iN)us as 
scattering data, where 3M G R is a sample of a random vector whose 
elements are independent Gaussian random variables with zero mean and 
unit standard deviation. 

Next, the approximation of Uoo{9'{^j),9{^j),k{^j)) is generated from 
Ua{xi,9j,kj) via numerical quadrature using (6). Note that the values 
Ua{xi,9j,kj) are only reliable for xi € D{0, 1). Nevertheless, this is not a 
problem in (6), since V is supported in D{0, 1). 

We have chosen R — 2.1 and as test potential V = Xi + 1-2X2 where 
Xi is the characteristic function of the annulus 0.7 < |(a;i,a; 2 )| < 1 and 
X 2 is the characteristic function of the square | 2 ;i| -|- |a; 2 | < 0.3 inside the 
annulus. 

5.2 The Born approximation: algorithm 1 

We first consider the fixed energy Born approximation. In this case, we 
can only compute 9{^) and 9'{^) for j^j < k/n. Combining the restriction 
1^1 < k/n and the fact that the mesh grid for ^ corresponds to the interval 
1^1 < N/{4R) we deduce that we have only far field data for the whole 
meshgrid in (19) when 



n 


If we choose k = IQ for example, this bound is 77 ~ 27. Thus, the far field 
data can be computed in the whole mesh with j £ when N < 27. 
The right hand side in (19) is assumed to be zero for > k/n. In Figure 
1 we show the real and imaginary parts of the reconstruction. Further, 
Figure 2 depicts the real and imaginary parts of this reconstruction, to¬ 
gether with the same for the true potential, in order that the reader can 
visualize the accuracy of the approximation. Only the values on the disc 
of radius 1.6 are shown. For comparison purposes, the same colour scale 
is used in all the pictures. 

Figure 3 shows the behavior of the error between the potential and 
the numerical approximation of the Born approximation in terms of k 
and N. Here, and in the sequel, this error is computed with the following 
approximation of the L^-norm, 



‘h 
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Figure 1: Numerical approximation of the Born approximation for fixed energy 
fc = 10 and N = 2^. The left figure corresponds to the real part and the right 
one to the imaginary part. 


real imag 


V 



VB,h 



Figure 2: Images of the true potential (top) and the numerical approximation 
of the Born approximation for fixed energy /c = 10 and N = 2^ (bottom). The 
left column corresponds to the real parts and the right column to the imaginary 
parts. For aesthetic reasons, only the values on the disc of radius 1.6 are shown. 
The colour scale is the same for the four pictures. 
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Figure 3: Numerical approximation of the Born approximation for different 
values of k. 


Note that this formula incorporates two types of errors: the one associ¬ 
ated to the numerical approximation of Vb and the distance of the Born 
approximation to the real potential V. There is no way to distinguish 
between these two errors since we do not know the Born approximation 
Vb for this particular example. 

From the results of Figure 3 we see that the error decreases for larger 
k and larger N. 

We now compare the rest of the Born approximations. 

In the case 2, fixed incident angle 8, we can only compute ^ outside the 
hyperplane ^-8 = 0. The right hand side of (19) is then assumed to be zero 
when ^ is in this hyperplane. The difference between the real potential 
and the Born approximation is very similar for any of the incident angles. 
Thus, we only show results for the case 0 = 7r/4 in Figure 4. 

In the case 4, full data, we have averaged the Born approximations of 
10 different incident angles equally distributed in [0,27r). In Figure 4 we 
compare the error between the potential and the real part of these Born 
approximations for different mesh sizes N x N = 2*'^ x 2^, M = 4, 5, 6, 7. 
We see that, for M > 6 the backscattering, fixed incident angle and full 
data generate the same error roughly. However, the fixed-A: case analyzed 
previously is significatively better than all these Born approximations. We 
also observe that, in contrast with the fixed-fc Born approximation, the 
error seems to stabilize after M — G. 

In the simulations above the synthetic far field data are computed 
in the same meshgrid as the one used to recover the potential. This is 
an important fact to avoid aliasing effects since we recover the potential 
from frequency data. If we use either a coarser or finer mesh to compute 
Uoo from the potential then it will contain necessarily an aliasing effect. 
To illustrate this, in Figure 5 we show the convergence of the numerical 
approximation for fixed k = 100 and backscattering, as N grows, when 
the data for the inverse problem (following Subsection 5.1) is computed in 
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Figure 4: Distance between the potential and the real part of the different Born 
approximation: d = 7r/4 fixed, backscattering and full data cases, depicted with 
circles, squares and crosses, respectively. 


a mesh twice finer than the mesh nsed to recover VB,h with the Algorithm 
1. We compare the error of the numerical Born approximation VB,h with 
respect to the real potential V in the mesh of size h, and with respect to 
a low pass hlter of the potential in the finer grid of size h/2 that removes 
half frequencies. As we see the latter error is smaller. 

This can be interpreted as the fact that we cannot recover anything 
better than a low pass filter of the potential V which takes into account 
the frequencies represented in the frequency mesh given by ^j. 

5.3 The iterative algorithm 

In this section we show the efficiency of the Algorithm 2 to approximate 
the potential. Note that in this case, at each iteration, we have to sim¬ 
ulate scattering data for the current approximate potential V" in order 
to recover the next approximation via the inverse Fourier trans¬ 

form. This is done following the idea explained in subsection 5.1 above 
for simulating the synthetic data in the experiments. 

In Figure 6 we show the difference between the real part of the ap¬ 
proximation and the potential for different values of fixed k and N = 2® 
grid points in each variable. We see that the algorithm decreases in the 
first few iterations and then it stabilizes. Thus, for hxed fc > 0 we observe 
a convergence V" —>■ F* yf F but close to F. 

We also observe in Figure 6 that, as k is larger, the approximation 
becomes better. The behavior is similar for other grid sizes. The con¬ 
vergence of the iterative method as fe —> oo is also illustrated in Figure 
7 where we Hx the number of iterations to n = 6 and compare F® — F 
for different values of k. We obtain numerically the following convergence 
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log N 



Figure 5: Distance between the numerical Born approximation of the potential 
and the potential when scattering data is computed in the same grid Vs^h — Vh, 
when a twice finer mesh is used to compute the scattering data Vb^h — 14 / 2 ) and 
when it is compared with a low pass filter of V in the finer mesh Vb^h — ^hl 2 ,LP- 
The left figure corresponds to the case k = 100-fixed while the backscattering 
case is in the right one. 


rate: 

error ~ (30) 

This rate is similar for other different initial data in . The situa¬ 
tion agrees with the estimate (18), obtained by Novikov in [18] for the 
modified sequence of approximated potentials IF”, in the sense that the 
convergence rate is a power of k. However, the exponent that we find 
numerically for F” — F is slightly better than the one given by estimate 
(18). In fact, Qoo = 0.5 when r = 4 while (30) holds for the I/°° potentials 
that we have tested. 

In Figure 8 we compare the convergence for the cases 2-4 with N = 2® 
grid points in each variable. We see that in the three cases the distance 
to the original potential decreases until a point from which it starts to 
increase. This is due to the fact that we are considering noisy scattering 
data and the algorithm converges to a “noisy” potential, which is not 
exactly the original potential V{x). 


6 Conclusions 

The inverse scattering problem for the Schrodinger equation consists of 
the recovery of the electrostatic potential F from scattering data measure¬ 
ments modelled by the far field pattern Uoo. The numerical approximation 
of this problem is investigated for different class of scattering data, namely, 
fixed energy, fixed incident angle, backscattering and full data. 

Two different algorithms in 2D are considered. The Algorithm 1 per¬ 
forms an inversion of the discretized far field pattern using the discrete 
Fourier transform. In this way the so-called Born approximation is ap¬ 
proximated. The Algorithm 2 goes beyond the Algorithm 1 through an 
iterative process solving successively the Lippmann-Schwinger (L-S) equa- 
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Figure 6: Distance between the potential and the real part of the numerical 
approximation when k is fixed, for different values of fc, in terms of the number 
of iterations. Here, the grid is taken with = 2^^ points. In the profiles, the 
circles, squares, crosses, rhombi and triangles refer to fc = 1, fc = 10, k = 100, 
k = 1000, k = 10000, respectively, as in Figure 3. 



Figure 7: Error in F" — V versus k after n = 6 iterations. The grid has 

jY2 _ 212 pQj]^^g_ 
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Figure 8: Distance between the potential and the real part of the numerical 
approximation in terms of the number of iterations, for fixed incident angle, 
backscattering and full data. An = 2^^-point mesh is considered. 


tion corresponding to the previous iteration result. The motivation of this 
iterative approach is the idea that these successive solutions to the L-S 
type equations approximate the scattered wave corresponding to the orig¬ 
inal potential, which is not proven and we propose as an open problem. 

A convergence result is provided for the Algorithm 1. This task is es¬ 
pecially difficult for the Algorithm 2 and is outside the scope of this work. 
The theoretical open problems that the Algorithm 2 pose, mentioned in 
the Introduction, are challenging and we hope they attract the attention 
of further articles. 

Both algorithms are tested from noisy scattering data and L°° poten¬ 
tials without rotational symmetry. 
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